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Abstract. Neural or cortical fields are continuous assemblies of mesoscopic models, also called neural masses, 
of neural populations that are fundamental in the modeling of macroscopic parts of the brain. Neural fields are 
described by nonlinear integro-differential equations. The solutions of these equations represent the state of activity 
of these populations when submitted to inputs from neighbouring brain areas. Understanding the properties of these 
solutions is essential in advancing our understanding of the brain. In this paper we study the dependency of the 
stationary solutions of the neural fields equations with respect to the stiffness of the nonlinearity and the contrast 
of the external inputs. This is done by using degree theory and bifurcation theory in the context of functional, in 
particular infinite dimensional, spaces. The joint use of these two theories allows us to make new detailed predictions 
about the global and local behaviours of the solutions. We also provide a generic finite dimensional approximation 
of these equations which allows us to study in great details two models. The first model is a neural mass model of 
a cortical hypercolumn of orientation sensitive neurons, the ring model |40|. The second model is a general neural 
field model where the spatial connectivity is described by heterogeneous Gaussian-like functions. 

1. Introduction. Neural or cortical fields are continuous assemblies of mesoscopic mod- 
els, also called neural masses, of neural populations that are essential in the modeling of 
macroscopic parts of the brain. 

They were first studied by Wilson and Cowan l43l [2) and are fundamental in the design 
of the models of the visual cortex proposed by Bressloff (10). Neural fields describe the mean 
activity of neural populations which are described by nonlinear integro-differential equations. 
The solutions of these equations represent the state of activity of these populations when 
submitted to inputs from neighbouring brain areas. Understanding the properties of these 
solutions is important for advancing our understanding of the brain. 

Among these solutions, the persistent states (or stationary solutions) are interesting for at 
least two reasons. First, when dealing with autonomous systems, looking for persistent states 
helps to understand the dynamics because they are an easy way to divide the phase space into 
smaller components. Moreover, as the connectivity is often chosen symmetric (see the review 
by CDO), the dynamics is described by heteroclinic orbits. Second, they are thought to model 
well the memory holding tasks on the time scale of the second which has been demonstrated 
by experimentalists on primates [13, 25l l36l . 

There are several items that arise in the mathematical description of these equations. The 
domain of integration (typically a piece of cortex) which is of dimension d and can be 
bounded or unbounded, the type of the nonlinearity (sigmoidal, Heaviside), the type of the 
connectivity function that appears in the integral (homogeneous, i.e. translation invariant, or 
heterogeneous), and the number p of neuronal populations that are modeled. 

The nonlinearity in the neural field equations is most of the time a Heaviside function 
which leads to several mathematical difficulties, one of them being being to specify the cor- 
rect functional space. We found easier in ||231 to use smooth nonlinear functions instead, to 
study the persistent states. Only a few papers use sigmoidal functions, e.g., ||3]|42], or an 
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alpha functiorQ |34|. In this paper we make the assumption that the nonlinearity is a sigmoid 
function, i.e. infinitely differentiable. Also, in this paper, the only assumption on the con- 
nectivity kernel is that it is square integrable, hence it can be heterogeneous (not translation 
invariant) as it is often assumed. 

Surprisingly, there are few papers dealing with the persistent states. Their authors use 
two main methods: Turing patterns (or bifurcation theory) and reduction to ODEs and PDEs. 
Bifurcation theory helps to understand the local structure of the set of persistent states and also 
gives the local dynamical structure. The first method is used in almost every paper, e.g. 11311131 
l5l [T0ll42l . using a translation invariant connectivity function (hence a convolution) and leads 
to the descrition of a lot of behaviours depending on relations involving the spatial frequency 
k and the time frequency u; traveling waves, breathers, persistent states. . . The numerical 
computation of these states depends on their stability. We generalize this approach as follows: 
we do not care about the spatial structure of the cortical states (indexed by k in previous 
studies) but consider a cortical state V(-) as a point of a (functional) vector space, and how 
this point varies (bifurcates) when the relevant parameters in the neural field equations vary. 
By doing so we are able to harvest a lot more results that do not depend upon the translation 
invariant assumption but of course also apply to this case. 

The second method is to reduce the persistent state equation to the problem of finding 
a homoclinic orbit to ODEs [33] when il = R and p = 1. If = W 2 , it reduces to the 
finding of homoclinic orbits to PDEs |34| . The obvious advantage is that they can use finite- 
dimensional tools 1 32 , 26 1 or PDEs methods, e.g. |3Q| . for the bifurcation analysis. Hence, 
the authors are able to compute the persistent states independently of their stability and to 
show their number as a function of the strength of the connections. 

Working with an unbounded domain is not biologically relevant and also raises some 
mathematical difficulties. We work with a bounded f2. Besides its biological relevance, 
this hypothesis is crucial for our mathematical analysis: it implies that there is at least one 
persistent state for any set of parameters and also it provides bounds for these persistent states. 

In a previous paper [23 1, we started the analysis of the neural fields equations defined 
over a finite part of the cortex from two different viewpoints, theoretical and numerical. We 
proved some results concerning the dynamics and gave a method to efficiently compute the 
persistent state (or stationary solution) under the assumption that it was unique. 

In the same article, we unsurprisingly found that the dynamics was poor (every initial 
condition converged to a single persistent state) if the stiffness of the nonlinearity was small, 
for example the system could not exhibit oscillatory behaviours. More importantly we did 
not give a way to compute the persistent states when the system featured more than one of 
such state. 

In this paper we relax the hypothesis on the uniqueness of the persistent state and try to 
understand the structure of the set, noted B, of stationary solutions as the parameter^] vary 
over large scales, hence we are not only interested in local behaviours. The local structure 
of the set B is understood using bifuraction theory (see 11351130 , 28 1) whereas degree theory 
helps to understand its global structure. This latter theory predicts stationary solutions that 
cannot be predicted using only bifurcation theory. Notice that bifurcation theory also gives 
the local dynamics. In order to compute numerically these persistent states, we propose 
a multiparameter continuation scheme that allows to compute non-connected branchesHof 
persistent states The case when there is no inputs can be (at least locally) done analytically 



1 e -1 /" H{v), H the Heaviside function. 

2 We choose the stiffness of the nonlinearity and the the contrast of the external inputs but our analysis applies 
to other parameters. 

3 A branch is a one-dimensional set of stationary solutions obtained by varying one parameter. 
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and is used to predict and verify the outcome of the numerical computations. 

For conducting our numerical experiments we consider a very general class of approx- 
imating connectivity kernels, the Pincherle-Goursat kernels ATI , which reduce exactly the 
dynamics to a system of ODEs whose dimension is directly related to the level of approxi- 
mation and can be arbitrarily large. Hence our choice to use infinite dimensional techniques 
to the integral equation is guided by the fact that it offers a simple, albeit abstract, concep- 
tual framework in which the behaviours of interest to us can be described in a clean and 
dimension-independent manner. When it comes to numerical experiments we use the system 
of ODEs provided by the Pincherle-Goursat kernels. 

The paper is organized as follows. In section [2] we introduce a very general functional 
framework for studying the neural field equations which allows us to pose the problem as 
a Cauchy problem in a functional space and to derive a number of useful properties for the 
(stationary) solutions of these equations. In section[3]we combine the results of the previous 
section with a bifurcation study in order to obtain more information about B. A numerical 
scheme is proposed to compute the structure of B. In section [4] we show how to reduce the 
neural field equations to a set of ODEs with an arbitrary precision through the use of the 
Pincherle-Goursat kernels. In section [5] we study in detail a neural mass model that reduces 
exactly to a finite set of ODEs, the ring model. In section |6| we compute the stationary 
solutions for a model of 2 populations of neurons on a bi-dimensional cortex connected by 
heterogenous Gaussian-like functions: it allows us to provide examples of the predictions 
obtained in sections |2] and [3] 

2. General framework. We consider the following formal neural field equation defined 
over a bounded piece of cortex and/or feature space £1 C R d . We wish to encompass in the 
same formalism the geometry of the cortex, seen as a bounded piece of R d , d = 1,2,3, and 
the geometry of such feature spaces as edge or motion orientations (directions) which are 
represented by a value between and it (2tt), or scale which is represented by a positive 
number, e.g. the spatial frequency. In section [5] we analyze in detail an example where the 
focus is on edge directions while in section|6]we analyze another example in which the focus 
is on the two-dimensional geometry of the cortex. 

V(r,t) = -L-V(r,t) + [3(t)-S(X(V(t)-9))}(r) + I ext (r,t) t>0 
V(-,0) = V (.) (ZA> 

This equation is an initial value problem that describes the time variation of the p-dimensional 
vector function V defined on il, starting from the initial condition Vq, a function defined on 
f2. At each time t > V belongs to some functional space, in effect a Hilbert space T, that 
we describe in the next section. We now discuss the various quantities that appear in ( |2.1[ ). 
J(t) is a linear operator from T to itself defined by: 

[J(i) • V(i)](r) = / J(r, r', i)V(r', t) dr', (2.2) 

where J(r, r', t) is a p x p matrix that describes the "strength" of the connections. We also 
describe in the next section the functional space to wich J(i) belongs and the conditions it 
must satisfy in order for the equation ( |2.1| i to be well-defined. 
The external current input, I ex t(-, t), is in T for all t > 0. 

The function S : R p -► R p is defined by S(x) = [S(xi), ■■■ , S(x p )} T , where S : R -> 
(0,1) is the normalized sigmoid function of equation 

S(z) = — 1— . (2.3) 



4 



R. VELTZ AND O. FAUGERAS 



It is infinitely differentiable on W and all its derivatives S^ 9 '(x), q = 1,2, • • • are bounded. 
For all integer q > 1 we note S^(x) the pxp diagonal matrix diag(S'^- ) (a;i), • • • , S^ q \x p )). 
Because of the form of the function S, the qth order derivative of S at x E M p is the multilinear 
function defined by 

D«S(s) ■(»!,••• ,y q ) = S^(x)- (yj ■■■y q ), w e R* i = 1, • • • , q (2.4) 

where y\ ■ ■ - y q is the component pointwise product of the q vectors yi, • • • , y ? of M p , i.e. 
the vector of MP whose fcth coordinate, k = 1 . • • ■ , p is equal to the product of the q fcth 
coordinates of each vector j/j, i = 1, • • • , g. 

A is the pxp diagonal matrix diag(A l7 • • • , X p ), Aj > 0, i = 1, • • • ,p that determines 
the slope of each of the p sigmoids at the origin. We note A m the maximum value of the AjS. 

9 is a p-dimensional vector that determines the threshold of each of the p sigmoids, i.e. 
the value of the membrane potential corresponding to 50% of the maximal activity. 

The diagonal pxp matrix L is equal to diag(^-, • • • , ^-), where the positive numbers 
Ti, i = 1, • • • ,p determine the exponential decrease dynamics of each neural population. 

As recalled in the introduction and detailed in 12T1 l22l this equation corresponds to a 
mesoscopic description of population of neurons which is called voltage-based by Ermentrout 
0181 . V in equation ( |2.1[ ) therefore has the biological interpretation of an average membrane 
potential of p populations of neurons. As pointed out by the same authors, there is also an 
activity-based mesoscopic description which leads to the following initial value problem: 

J A(r,t) = -L tt .A(r,f) + S(A([J(0-A](r,t)+I rat (r,f))) i>0 

\ A(,0) = Ao(-) ^ } 

L a ^ L (see ifTHl ) because they do not have the same biologiocal meaning: One is related 
to the synaptic time constant and the other to the cell membrane time constant. We let L Q = 

diag(«i, • • • ,Op). 

The two problems ( |2.1| i and ( |2.5| l are closely related. In particular there is a one to one 
correspondence between their equilibria, as recalled below for the ring model of section |5j 
which is an activity-based model. 

2.1. Choice of the appropriate functional space. The problem at hand is to find an 
appropriate mathematical setting, i.e. to choose the functional space T, for the neural field 
equations based on three criteria 1) problems ( |2.1| i and ( |2.5| l should be well-posed, 2) its 
biological relevance and 3) its ability to allow numerical computations. 

A Hilbert space is appealing because of its metric structure induced by its inner product. 
Hence, a natural choice arising from the definition of the linear operator J would be the space 
h 2 (il,W). But this space allows the average membrane potential to be singular which is 
biologically excluded. For example the function r — ► | r | 1 / 2 g L 2 (fl, MP) if d > 0. It would 
be desirable that this potential be bounded on the cortex. A way to achieve this is to allow 
for more spatial regularity of the membrane potential: it is reasonable, E.g., from optical 
imaging measurements, to choose r — > V(r) as being differentiable almost everywhere. 

For technical reasons that will become clear later we choose the Sobolev space 

t = w m ' 2 (n,M p ), men, 

with the inner product: 

m 

(Xi,X 2 ) y = (^ Qx i^ QX 2>l W p) > VXx, X 2 e T, (2.6) 

\a\=0 
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where, as usual, the multi-index a is a sequence a = (ot\, ■ ■ ■ , ad) of positive integers, 
\a\ = a i-> anc ^ tne symbol D a represents a partial derivative: 

^QiH Va d 

° atp= dr?-.9r?' P > 

where ip : O — > M and r = (n, • • • , r^). 

Note that .F = VF m ' 2 (0) x • ■ ■ x W m ' 2 (Q). Later we use the notation Q for VF m ' 2 (0), 

s J 

V 

i.e. JF = gP. 

2.1.1. Choosing the value of m. The value of the integer m determines the regularity of 
the functions V that represent the cortical states as well as that of the connectivity function J. 
It turns out that, depending upon the relative values of m and the dimension d — 1, 2, 3 of the 
space where we represent the cortical patch il and/or the feature space, T is a (commutative) 
Banach algebra for the point wise multiplication. This is important when using bifurcation 
theory because we need to use Taylor expansions of the right-hand-side of ( |2.1[ ). 

Being a Banach algebra property requires some regularity of f2, i.e. of its boundary 
<9fL Because in the numerical experiments of section|6]we work with the open square O = 
[—1, l] 2 , we will assume the weakest regularity, the so-called cone property, see U Chapter 
IV] for a definition. In practice one can safely assume that d£l is C 1 . Roughly speaking it 
means that this boundary is a C 1 -manifold of dimension d — 1 of M. d . The exact definition, 
called the uniform C 1 -regularity, can again be found in [1, Chapter IV]. It implies some 
weaker regularity assumptions, such as the cone property. 

Under this assumption on il we can adapt the theorem in |fl] Chapter V, Theorem 5.24] 
to obtain the following result : 

COROLLARY 2.1. If the relation 2m > d holds, then T is a commutative Banach 
algebra with respect to component pointwise multiplication, i.e. for all U and V elements 
of T their component pointwise product UV is in T and there exists a positive constant K* 
such that 

||UV||^ < K* HUII^ ||V||^ VU,Vef, 

where UV = (U x Vi, • • • ,U p V p ). 

Let us explore some consequences of this proposition on our possible choices for the 
value of m. These choices are guided by the idea of constraining the functional space T as 
little as possible, given the fact that W m ' 2 C W"' 2 for all integers < n < m 
d = 1: The relation 2m > 1 implies m > 1 and we choose m = 1. 
d = 2,3: The relation 2m > d implies m > 2 and we choose m — 2. 
Similar results hold for higher values of the dimension d but these cover the cases discussed 
in this article. 

We summarize all this in the following proposition 

PROPOSITION 2.2 (Choice of J 7 ). If O is has the cone property, in particular if it is 
uniformly ^-regular, then we have T = W 1 ' 2 (0, W)for d=\andT ' = W 2 < 2 (0, MP) for 
d = 2, 3. In all three cases, T is a commutative Banach algebra with respect to component 
pointwise multiplication. 

2.1.2. The choice of J. We assume that J(-, -A) e W m > 2 (0 x fl,RP x P) for all t > 0. 
As a consequence, the Frobenius JF-norm of the linear operator W(t), noted || 3(t) \\ ^, is 

well-defined 



2 V- ' 2 



H,M=o 



||J(t)|£= X / D«D a J(r,r>,t) 
Inxn 



dr dr', 

F 
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where it is understood that the partial derivative operator D a (respectively D a ) acts on the 
variable r (respectively r') and || J(r, r', t)\\ F is the Frobenius norm of the matrix J(r, r', t), 



\\3(r,r',t)\\ 2 F = J2M^r',ty 



This hypothesis also ensures the existence of the operator J as a linear operator from T 
to T. 

PROPOSITION 2.3. Assume that J (-,-,£) S W m > 2 (il x (l,W x P) for all t > 0. Then 
equation ( |2.2| l defines a linear operator from T to T. 
Proof, see appendix |A|D 

Since S is a W -valued bounded function on fi and such functions are in L 2 (f2) the right-hand 
side of ( |2.1| i is an element of T. Similarly, because S is also infinitely differentiable and all 
its derivatives bounded, the right-hand side of ( |2.5| l is an element of T. 
Finally we note || | J(i)| || the operator norm of J(i), i.e. 



sup IIV 



It is known, see e.g. [31 ], that 

ll|J(*)lll<l|J(*)IU 

Remark 1. Note that the relation 2m > d and our regularity assumption on f2, imply, 
through the so-called Sobolev imbedding theorem, [1. Chapter V, lemma 5.15], that 

where Cg(fl) is the set of continuous bounded functions on f2. The consequence is that the 
state vectors V and the connectivity matrix J are essentially bounded, which has certainly 
the right biological flavor. 

2.2. The Cauchy problem. In this section we prove that equation ( |2.1) is well-posed 
and provide some properties of its solutions. 

We rewrite it as a Cauchy Problem, i.e. as an ordinary differential equation on the func- 
tional space T. This turns out to be convenient for the upcoming computations. 

d J = -L-V + R(t,V) t>0 
V(0) = v e^ ; 

The nonlinear operator R is defined by 

R(t, V) = J(t) • S(A(V - 8)) + l ext {t) (2.8) 



Proposition 2.3 shows that R(t, •) : T — > T for all t > 0. We have the further properties: 

LEMMA 2.4. R satisfies the following properties : 

• V integer q>0, R(i, ■) G C q (T,T) and D q R(t,V ) = A 9 J(t)S^^(A(V — 0)) 
for all Vo in T. 

• '||R(t,Ui) -R(t,U 2 )||jr < A m ||J(t)||j7 1| Ui - \5 2 \\jr for all t > and for all 
Ui, U2 in T. 

• R(t, •) is a compact operator on J- for all t > 0. 
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Proof. It is easy to see from the definition ( |2.8| ) of R that, if it exists, D q H(t 1 Vo) [Ui , • • • , U g ] 
X q 3(t) ■ (S^(A(V Q - 9)) • (Ui • • • U,)). The notation U x • • • U g is the same as in the def- 
inition of D q S in equation Q2.4) , i.e. the component pointwise product of the q functions 
U lr -- ,V q ofF. 

The g-multilinear operator D q R(t, V ) is well-defined because, according to corollary 



2.1 Ui • • • XJ q is in T . The first property follows immediately. 

It remains to show that the g-multilinear operator D q H(t, Vq) is continuous. We apply 



corollary 2.1 once more to show that 



\\D q R(t,V ) ■ (Ui • ■ ■ U q )\\ F < ||A 9 J(i)S^(A(V - 0))||^ ||XJ X • ••UJ^ < X J] ||U 4 

i 

for some positive constant K, that is D q H is continuous and R(t, •) G C q {T , T). 
The second property is proved in l23l . 

For all a, | a | < m U — > 9" J ■ U are compact operators on L 2 (fi, MP) because they are 
Hilbert-Schmidt operators (see B4l chapter X.2]). Hence for any bounded V n sequence in 
T hence in L 2 (ft, W), there exists a subsequence such that D a 3(t) ■ S(A(V^ - 0)) are 



convergent in L 2 (f2, W) for \a\ < m. Then R(t, V^) is convergent in T because of (2.6 1. 
We have proved that R(f , ■) is a compact operator on T . □ 

Following closely EH we have the following proposition. 
PROPOSITION 2.5. If the following two hypotheses are satisfied: 
1. The connectivity function 3 is in C(K + ; W m ' 2 (f2 x fl,R pxp ) for the values of m 



given in proposition 2.2 and bounded, HJ^)!)^ < J, t > 0, 
2. the external current I cxt is in C(M. + ; T), 
then for any function Vo in T there is a unique solution V, defined on R + and continuously 
differentiable, of the initial value problem \2.\\ . This solution depends upon 'dp parameters, 
the slopes A, the thresholds 6 and the diagonal matrix L. 

Even if we have made progress in the formulation of the Neural Field equations, there 
still remains the unsatisfactory possibility that the membrane potential becomes unbounded 
as t — > oo. However this is not the case as shown in the next proposition: 

PROPOSITION 2.6. If the external current is bounded in time ||I C xt||jr < / C xt> for all 
t > 0, then the solution of equation ( |2.7| > is bounded for each initial condition Vo G T. 

Proof. Let us define / : E x T -> E+ as 

f(t, V) d ^ f (-L • V + J(i) • S(A(V - 9)) + I ext (t),V) F = 1 d||V "^ 



lr 2 dt 



We note T max — maxj=i... , p Tj and notice that 



/(i,V)<--^||V|| 2 r + (J + / erf )||V||^. 



Thus, if ||V||^ > 2r max (W + I ext ) d = R, f(t,V) < -2r max (W + I ext ) 2 = ~5 < 0. 
Let us show that the open ball of T of center and radius R, Br, is stable under the 



dynamics of equation (2.1 1. We know that V(t) is defined for all t > 0s and that / < on 
dBji, the boundary of Br. We consider three cases for the initial condition V . 

If V E B R and set t = sup {t \ Vs e [0, t] , V(s) € Br}. Suppose that r G K, then 
V(r) is defined and belongs to _B#, the closure of Br, because Br is closed, in effect to 
dB R . We also have ^||V|||-| t=T = /(r, V(r)) < -J < because V(r) G Thus we 
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deduce that for e > and small enough, V(r + e) G _B/? which contradicts the definition of 
r. Thus r ^ R and _B/j is stable. 

Because / < on dB R , V G <9-Br implies that Vt > 0, V(t) G B R . 

Finally we consider the case Vo G £>B R . Suppose that Vt > 0, V(i) £ B R , then 
Vi > 0, Jj||V||^r < —28, thus ||V(i)||jr is monotonically decreasing and reaches the value 
of R in finite time when V(i) reaches dB R . This contradicts our assumption. Thus 3t > 
0\V(T)eB R .U 

Corollary 2.7. 7/V £ B R andT = inf {t > 0|V(i) £ 77?en 

r< l|Voll^-fi 2 
26 



This proposition shows that B R is an attracting set and that it suffices to study the dy- 
namics within this set to have the long time behavior of the solutions of the Neural Fields 
Equations. 



V{, of (|2.1 



This attracting set contains the stationary solutions of (2.1 1, we devote the next section 
to their study. We quote a result from ll23l concerning their stability: 
Proposition 2.8. If the condition 

X m p(3 s ) < 1, 

holds, then every stationary solution of \2.1^ is globally asymptotically stable. X m — max A,, 
J s is the symmetric part (J+J*)/2 of the operator J, and p(J s ) its spectral radius. We define 
X L to be p(3 s )~ 1 . 

Similar results hold for the activity-based model \2.5\ . 

2.3. Pro perties of the set of persistent states . We look at the equilibrium states, noted 
when I ext and J do not depend upon the time. Our goal is to estimate their 
number and, if possible, to compute them numerically, for a given set of parameters. 

It is quite demanding to do it at a given point in the parameter space except in some very 
special cases^] We note that when A = (or J = 0), the stationary equation is trivially solved. 
Hence, we can think of deforming this trivial solution to a solution when A ^ 0,J ^ 0. This 
raises a number of questions. Does such a "manifold" of solutions exist i.e. can we link the 
trivial solution to a solution for any given set of parameters? If yes, are there any other solu- 
tions? How do these "manifolds" look globally? These questions concern global properties of 
the set of solutions (existence of branches, existence of intersection points, connectedness. . . ) 
are difficult to answer. We provide some partial answers in the remaining of the article. 

Before going deeper in the analysis, we need to simplify the parameter space. The equi- 
libria, also called bumps, or persistent states , are stationary solutions (independent of time) 
of 

= -L ■ v{ + J ■ S(A(V{ - 6)) + I ext , 
We redefine J as L _1 J, V as V — 9 and l ext as L _1 ■ I ext + and restrict our study to: 

0= -V{+J-S(AV{)+I cxt (2.9) 



Still, equation (2.9 1 contains many parameters such as the ones describing J and I cxt , or 
the slopes A. Which parameters to choose for the continuation method: A or J ? We decide 
to set J and control A for two reasons : 



4 When the nonlinearity J ■ S(A(V^ — 8)) is small compared to the linear part L • V^, we know there exists a 
unique solution and how to compute it efficiently. This was proved in [23 j . 
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• the stationary solutions are bounded for A £ R5_, see proposition 



2.9 



1 which is not 



»■+. 

the case when || J|| — > oo. This proves to be useful numerically. 
• previous studies usually use a Heaviside nonlinearity which is formally equivalent 
to our nonlinearity when A'='oo, varying A can thus bridge the gap with previous 
approaches. 

As a matter of fact the techniques we are about to expose are applicable to any set of param- 
eters with minor modifications. Hence, we now focus on the influence of the slopes A on the 
solutions of (2.1 1. We make the assumption that they are all are equal to A, A = AId p , A > 0. 
where Id p is the p x p identity matrix. The equation becomes 

= -V{ + J • S(AV{) + I ext d ^ f -F(V{, A) (2.10) 
It is clear that when A = 0, the stationary equation is trivially solved by 

V f Q = J-S(0)+I ext = ^J-l + Iext, 

where 1 is the p-dimensional vector with all coordinates equal to 1. Let B\ be the set of 
solutions of equation ( |2.10| i for a given slope parameter A : 

S A = {V|F(V,A)=0} 

We next provide some properties of the sets B\. 
Proposition 2.9. 



1. The persistent states satisfy the following inequality 

rf 



V{ - V} 



< II J I 



cxtlljr 



where So : M — » M is the "shifted" sigmoid defined by Sq(x) = S(x) — S(0) and 
the constant B\ is defined in proposition \D.2\of appendix\D. 1\ 



2. If the condition 



A||J||^< 1 



(2.11) 



is satisfied, then #£>a = 1. We define A* to be ||J||^- . 

3. VA e R+ B x ^ 0, 

4. If we know an even number of solutions in B\, then B\ contains at least one more 
solution. 

5. Let < a < b be two reals, and consider the set B = U^efa, b]B\ x {A}. Then B 
contains a connected component C which intersects B a x {a} and Bf, x {&}. 

Proof. 



1. From lemma 
Therefore 



D.l 



in appendix|D. 1 |we have S (XV^) 2 < S {X 2 (V^) 2 ), i = 1,- • • ,p 



'A2\ ■ _ 



A 



2 P 



i=l i=l \ P i=l J 

The second inequality comes from Jensen's and the fact that Sq(-) is concave in R + . 
It then follows, using again Jensen's inequality and the fact that So is monotonously 
increasing, that 



L 2 (f2,Rp) 



<p\n\s 



A 2 

P \n\ 



v 



l 2 (o,: 



<p\tt\S 



P \n\ t, 
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Now 



A v 



J-S (AV{) = p Q J-S (AV 



< 



S (AV{) 



L 2 (fi,] 



L 2 (n,Rp) 



< 



S (AV{) 



l«|=o 



L 2 (f2,HP) 
2 

L 2 (n,Rp) 



(2.12) 



The inequality then follows from propositio n|D.2| 

Use the Picard Theorem. As shown in figure |2l2|this imposes that A* < Aj> Indeed, 



-l 



and HIJ.JI = II I Jl II, and II I J| II < ||J| 



we have A* = J 



T 



< 



as p(3 s ) < 
X L = p(3 s ) 

The first property is that it is non empty: in ||231 we proved that persistent always 
exist in L 2 (f2,R p ) for any positive values of A. However, our ambient functional 
space is different since we require more space regularity. Let us consider a per- 

But 



R(v{). 



sistent state V{ £ L 2 (fi,R p ) given by ||23]|. It satisfies V 
R : L 2 (f2,]R p ) — ► T because of the assumptions on J(r,r') and I ex t- Hence any 
persistent state in L 2 (f2. R p ) in fact belongs to T. 

Suppose that B\ has an infinite number of solutions, then the proposition holds. If 
now, B\ has a finite number of solutions, we can assume that these points are non 
critical. Then according to the Leray-Schauder degree theory sketched in appendix 
|B]we have 

deg LS (^(V, X),B r , 0) = £ si § n dct LS ( J D v F(V{, X),B r , 0) = 



v{ee A 



£ signdet LS (Id-AJ^S(AV{)), 



where r = 2B\ and B\ is defined in proposition [1X2] in appendix |D.1| We prove in 
corollary B.2 in appendix [B] that the first term is equal to 1. Suppose now that B\ 
contains an even number of points, say 2k among which I correspond to a negative 
sign and hence 2k — I to a positive sign. The sum that appears in the last term is 
equal to 2k — 21, hence even. Hence B\ must possess an odd number of points. 
5. The proof uses the Leray-Schauder theorem, see appendix [B] We apply the theorem 
to the function F : T x J — > T which is of the form Id + to, with to(-) = — R. 
Because to is compact on T x J (see proof in [21 1), B a x {a} bounded, and there ex- 
ists an open bounded neighbourhood U a of B a such that deg LS (F(-, a),U ai 0) ^ 

the conclusion follows since the connected compo- 



(corollary B.2 in appendix |B| 



nent cannot be unbounded, B being bounded. 



This proposition answers some of the previous questions. For any positive value of A, 
there is always at least one persistent state and we can find a way to connect the trivial solution 
Vq to a persistent state corresponding to an arbitrary value of the parameter A. We return to 
this connection later. All this is true if we choose, say || J||jf, as a parameter instead of A. We 
will also see that not all the solutions in B\ are in the connected component of (Vq,0) in 
T x M+. 
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Regarding the connection between Vq and V, proposition 2.9 does not give us any 
indication on its regularity but we have the following corollary. 

COROLLARY 2.10. Let a and b be as in proposition \2.9\ For all e > there exists 



V,- V,: 



-ill^ 



< e for 



a finite sequence (Vx,Ai), (V n ,A„) of points of C such that 
i = 1, • • • , n — 1 and Ai = a, A n = b. 

Proof. C is connected for any topology equivalent to the product topology of T x [a, b], 
e.g. for the metric defined by d((V l7 Ai), (V 2 , A 2 )) = ||V"i - V 2 11^+ |Ai - A 2 |. Since it is 
connected for this metric, it is also well-chained 1121 . and the conclusion follows. □ 



In fact, except at points where the Jacobian of F is non-invertible (such points like 
B in figure |2.1| are potential bifurcation points), the impl icit functions theorem tells that 

5 imposes strong constraints 



2.9 



A — > (V^, A) is differentiable. Hence in effect proposition 
on the set B\ as shown in figure 2.1 The horizontal axis represents the parameter A, the 



vertical axis the space T where the solutions of (2.10i live. The curves represent possible 
solutions as functions of A. The configuration in the lefthand part of the figure are forbidden 



by proposition 2.9 5 while those on the righthand side are allowed, the green curve being an 
example of a continuous curve s — > (Vw g N, A(s)) from [0, 1] to T x [a, b]. 



solutions of F(V,A 




FIG. 2.1. In the lefthand part of the figure, there is no connected curve of solutions in [a,b\: this is forbidden 
by proposition \2.9\ 5 which states that we must be in the situation shown in the righthand part of the figure, where 
the green curve connects B a X {a} and Bt X {b}, see text. 



Proposition 2.9 5 gives a very interesting general (non-local) property of the set of solu- 
tions. But it is non-constructive, for example it does not tell us which branch to chose at point 
B in figure 2.1 and we need to compute all the branches to know the path to A = b. Hopefully 
such branching points as B are very rare and one only sees turning points rather than branch- 
ing points (such as the one at Ai in figure 2.2 left): any perturbation will indeed destroy such 
branching points (see figure 2.2 right). Hence, if one continues the trivial solution (obtained 



for A = 0), one will typically find a curve like the green one in figure 2.2 Right. This may 
lead to the wrong conclusion that for A big enough there is only one stationary solution in- 
stead of three. The problem is to find a way to compute, if it exists, the second, red, curve 
which is not connected to the green curve, hence not attainable by A-continuation. 

An idea, directly suggested by the above picture is to restore the branching points by 
perturbation. Among all possible perturbations, we choose one of the simplest, we vary the 
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contrast el ex t by varying e from to 1, as suggested by experimentalists who usually provide 
neural responses as functions of the contrast. 

The conclusion is that if we want one solution for A ^ 0, we use a A-continuation of 
the trivial solution, but if we want more than one (or the maximum number of) solutions for 
A 0, then we perform a (A, e) -continuation of the trivial solution. 




FIG. 2.2. The white zone is the domain of existence o/VA The grey zone is excluded thanks to proposition 
\2.9\ 1. Left: non generic situation, a transcritical bifurcation occurs at X\. Right: same as left with a small pertur- 
bation (like a tiny change of the external input I cx tj, the transcritical bifurcation has opened up. \l is defined in 
prop \2.8\ and A* in prop \2~9\ 

Remark 2. An interesting question is to predict how close to \l can be the smallest value 
of X where a "turning point" occurs. 

However, when performing this (A, e) -continuation of the trivial solution, we may end- 
up with a lot of data. Hence, we need to have an idea of of the structure of the set (V{ , A, e) . 
The idea is, again, to think of this set as a deformation of an easier to compute set of stationary 
solutions. This is done in the next section. 

Remark 3. It is highly possible that this (A, e) -continuation scheme still misses some 
solutions. One possibility is to introduce a third parameter, but it may become quickly numer- 
ically untractable. 

3. Exploring the set of persistent states. We exploit the fact that in the neural field 
equation ( |2.10| > the ratio between the external current I ext and J is not fixed a priori. Hence 
when studying the NFE, one would rather look at 

- V + J-S(AV) + eI 0Xt = (3.1) 

where e > 0. The persistent states now depend upon the pair (A, e) i.e. V{ . The idea is to 
infer the persistent states of from those, Va.o of 

-V + J-S(AV) = (3.2) 

We further simplify the problem by considering 

- V + J-S o (AV) + /iJ-S(0) =0, (3.3) 



where So is defined in proposition 2.9 We recover equation ( |3.2| i when p, = 1, The advantage 
is that when p = 0, we can say a great deal about the persistent states. 
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3.1. A simpler case. This simpler case reduces to the study of the previous equation 
when fj, — 0: 



- V + J S (AV) = 



(3.4) 



In other words, we infer the persistent states of ( |3.1[ ) from those, Va,£=o.^=o of equation 
( |3.4| ). This case has been studied a lot by several authors lfl8l [141 [6] when a constant current 
I ex t is applied, which amounts to changing the threshold 9 in S. 

V = is a trivial persistent state for (|3.4|i. Recall that a necessary condition for the equa- 



tion F(V, A) = to bifurcate at the solution (V{, A) is that DyF(V{, A) is non-invertible. 
For technical reasons it is interesting to slightly modify the nonlinear term R in \2.1\ by 



subtracting from it its linear part as follows. We rewrite ( 3.4 1 as 



withR(A,V) = 0(||V||§r) and 



L A • V + R(A,V) = 



-Id + AJL»S o (0) = —Id + - J 



(3.5) 



(3.6) 



The operator La satisfies the following properties. 
Proposition 3.1. 

• The operator J is a compact linear operator on T. 

• L\ = —Id + j J is a Fredholm operator of index ( useful for static bifurcation, see 

mi 

• La is a sectorial operator (useful for dynamical bifurcation, see [28]). 



Proof. The first property was proved in lemma 2.4 



Because J is a compact operator, the kernel of —Id + h J is of finite dimension equal to 
the codimension of its image, hence its index is 0. 

La is a sectorial operator because it is a bounded linear operator on a Banach space, see 
e.g. (29]. □ 



We can therefore state that the values of the parameter A that determine the possible 
bifurcation points are: 



Ar, 



4 

0~n 



1, 2, 



where a n is an eigenvalue of the compact operator J. We assume in the sequel that Ai < 
A 2 < • • • . 

Additional properties have to be met in order to have a bifurcation point (see |35|) but 
they are satisfied if every eigenvalue of J is simple. We will assume that the first k > 1 
(k is arbitrary) eigenvalues are simple because this class of operators is dense in the set of 
compact operators set (see ll35l ). We denote by e n (respectively by e* ) the eigenvector of J 
(respectively of the ajoint operator J*) for the eigenvalue er„. 



A simple adaptation of lemma 2.4 shows that R is C 9 , for all integer q and we can 
consider (if it exists the minimal integer q > 2 such that 



0^4") = (Gg(e„,A n ),e*) 



According to lemma 2.4 



1 



G q (e n ,X n ) = -L> 9 R(A„,0)[e 



A« 



^ 9) (0)J-(e„---e„) 
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By parity, Sq (0) = s q ^ if and only if q is odd. Hence, the parity of q, which is odd, 
tells that we have a pichfork bifurcation at A„. In particular, the bifurcated branch is given by 



V{ = x(X)e n + o(x) 

Thus the bifurcation portrait is a set of branches C n emanating at points (at least for 



n < k) (0, A„). Thanks to proposition 3.1 we can apply the Lyapunov- Schmidt procedure 



(see |28"1 ) to get the following reduced equation oni£ 

= (-1 + A(j„/4)x + g(x, A) = (-1 + Xa n /4)x + X ^x q + o{x q ) 

Depending on the sign of Xq > tne pichfork branch is oriented toward A < A„ (resp. 
A > A„) if Xq 1 ^ > (resp. Xq < 0)- Let us note e q n the vector e n ■ ■ ■ e n . We have 

gtimcs 



(n) _ K s g ;j j .1 _ ^n s q I q j* . * \ _ K f? / q , , 

<> q. 



From a practical standpoint, the last inner product is difficult to compute since it requires the 
computation of the eigenvectors of the adjoint operator J* of J for the inner-product of the 
Hilbert space T = W m ' 2 which are in general different from those of the adjoint operator 
J2 of J for the inner-product of the Hilbert space W 0,2 = L 2 . This computation is greatly 
simplified by the following proposition. 

PROPOSITION 3.2. Let e* L2 be an eigenvector of the adjoint operator 3* L2 of 3 for the 
inner-product of the Hilbert space L 2 , associated to the eigenvalue A (assumed to be simple). 
Let e*-p be the corresponding eigenvector of the adjoint operator J^r of 'J for the inner-product 
of the Hilbert space T = W" 1,2 . Then it is possible to choose ej- such that the following 
holds 

(U,e»^HU,e! 2 ) 2 VUef 

Proof, see appendix|E] □ 

(n) 

In effect, in order to obtain the sign of Xq , we only have to compute the eigenvectors 
of 3 2 and to compute inner products in L 2 . 

Thus, we have found local branches of stationary solutions and continue them globally 
in order to obtain the global branches named C„ . An interesting question, yet unsolved, is 
to know whether the branches C n are connected. Some results exist in that direction in the 
line of those of Rabinowitz (see [38, 35]) but dont give much insights in the general case 
(d > l,p > 1...). However, they can be used to derive some properies of C\. 

Proposition 3.3 (Turning point property). Ifxq > 0. tnen 3At < Ai such that 



VA G (At, Ai), (3.4 \ has at least 3 solutions and At = min {A/ (V^, A) € Ci}. 

Proof. Let C\ be the connected component in B where B = {(V, A) | V 7^ 0, (V, A) satisfies ( |3.4| | 
to which (0, Ai) belongs. Then (see ||35l ) C\ is unbounded in T x M + or contains a point 
(0, A n ), n > 1. In either case, C\ exists until A = A2. 

C = C\ n {J- x [A*, A2]) is closed and bounded (because every solution V{ is bounded 
in J 7 ). Let us show that C is compact in T x [0, A2]. Consider a sequence (V n , s n ) in C. As 
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s n is bounded, we can assume it is convergent. We also have V„ = J • S(s„V„) + I cxt . As 
(V, A) — > J ■ S(AV) + I is a compact operator, there exists a subsequence (V^,( n ), SM n \) 
such that J ■ S(s0/ n )V^,( n \) + I is convergent, hence ~VM n ) lS converging. We have proved 
that C is compact. Hence LT K +(C) is a compact subset of K + . Then inf LT R +(C) is a rain 
written At € n R + (C). As it is an inf, there exists a sequence s n associated to a V n in C such 
that s n — > At- But as C is compact, we can assume that V„ — > Vt- Then (Vt, At) € C is 
called a turning poin^ 

So we have proved that At < Ai. But in the case q odd with > 0, Ci exists for 
A < Ai and At < Ai. 

Now, VA e (A T , Ai), 3(V{, A) € C x \ {(0, A)} because d intersects {0} x R + only at 



bifurcation points (0, A„) located 'after' Ai. Then because of proposition 2.9 part 4, there are 
at least three solutions. □ 

Remark 4. If the sigmoidal function had satisfied S^(0) ^ 0, we would have seen 
transcritical bifurcations and the previous proposition still holds in that case. 

Remark 5. If we were able to prove that the branches do not intersect, the previous 
proposition would apply to all branches C n satisfying the required conditions. 

Even if we do not deal with the dynamics, we can say a little using [28 , 35 , 32 1. At every 
bifurcation point, the center manifold at V = (see [28]) is one-dimensional while the di- 
mension of the unstable manifold increases as A crosses values corresponding to transcritical 
points. Hence for large As, we can have a large unstable manifold. Note that every value of A 
is biologically plausible because the locally (around V = 0) exponentially divergent dynamic 
is bounded (see proposition |2.6| l, which can make the "global" dynamics very intricate. In 
effect, when A grows to infinity, the sigmoid tends to a Heaviside function which acts as a 
threshold. 

Remark 6. Note that it is easy to characterize Breathers (see H37\\16\\24\l ) in this frame- 
work. It is sufficient to choose a connectivity function 3 with complex eigenvalues. For 
example with a model of two populations (p = 2), one excitatory, one inhibitory, such that: 



J(r,r') 



Gi(r,r') -G 2 (r,r') 
G 2 (r,r') Gi(r,r') 



Gi symetric > 0, [G±, G 2 ] = 



Thena(J) = {gi_ n ±ig2,m9i,n € cr(Gi)} and a Hopf bifurcation occurs at Xn.n = 2S'(a)g 1 
if the eigenvalue is simple. Because R is G 3 , we can also compute its first Lyapunov coeffi- 
cient (see [28],[32]p.l07) and find 

A 3 



, A H,n b 2 I 2 \ 2S 3 S! 

with u n = V 1; i + „ g '- , s 2 = 5( 2 ) (0), S3 = S (3) (0) 

3.1.1. The case d = l,p = 1. We can use the results of Rabinowitz (see 11381 ) in the 
case d = p = 1 and when the connectivity operator J is a symmetric convolution. In that case 
the zeros of the eigenstates of J are simple on £1 C K (because they are cos and sin functions). 
Moreover the Taylor formula allows to write Sq(x) = xh(x) with h > 0. Then from [38 1, 
it follows that the bifurcated branche^]C n (pichfork or transcitical) meet (A„, 0) and oo and 
are caracterized by the number of zeros of their elements ; namely V (V-^ , A) € C n , has 
exactly n — 1 zeros in f2. As a consequence, the branches C n do not intersect. 

5 Not in the sense of [32], here it denotes a local extremum in the parameter along the curve of solutions. 
6 Here they are connected components in C 1 (SI, IR) X M.+ rather than in T X BS+ . 
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If we were able to prove that all the stationary solutions are connected to the zero solu- 
tion, we would have completely caracterized B. Nevertheless, it is still a lot of information. 

Remark 7. The study of the simpler case is very important for the numerical computa- 
tions. Indeed, it allows analytical predictions. When one performs the (A, p, e) -continuation, 
one should look at the section of solutions (A, p = 0, e = 0) and compare to the predictions 
of the simpler case to see if the numerical did not miss some solutions ( that may happen when 
the system has some symmetries). 

3.2. Returning to the original equation. The overall picture that emerges from the 
previous section is interesting despite the fact that some of its features are hard to justify from 
the biological viewpoint: no external input, rate function So possibly negative. The reason is 
that it gives clues about the set of persistent states when e ^ an d provides a way to compute 
them numerically. Starting with the trivial solution Vq of (3.2 i when the slope parameter A 
is null, we can perform a numerical continuation with respect to the two parameters A, e. 

When e 7^ 0, the only bifurcations that are possibly unaltered are the turning points. The 
transcritical/pichfork bifurcations will be "opened" as described below. We are still able to 
predict the stability near the points (0, A„). Let us note / = (/iJ • S(0) + £l ex t, Ki) r- Then, 
the Lyapunov-Schmidt method |28 35 1 leads to 







-1 



Ao>, 



X^x" + I + o(x q ) 



(3.7) 



Notice that when q = 2, we find the same equilibria as those of the Bogdanov-Takens 
normal form and for q = 3, we obtain a cusp. Solving the polynomial equation ( |3.7| i allows 
us to describe different opening scenarios depending on the sign of I. The cases / > are 



shown on figure 3.1 The case of Ai is a little bit special according to proposition 3.3 and is 



shown in the righthand part of figure 3. 1 



00 A 




J<0 




FIG. 3.1. Opening of the pichfork bifurcation at the first (right) and nth (n > 1, left) eigenvalue. Black dots 
indicate saddle nodes. Note that there are three such points, see proposition \3.3\ for the first eigenvalue. Note that 
the branch may have a more intricate shape than the one shown. 



Remark 8. When the first eigenvalue generates a 'well oriented' pichfork branch, propo- 
sition 3.3 says that a turning point must occurs on this branch: there are two turning points on 



this branch. Still we have a local description and it would be interesting to have more global 
results for example concerning the behaviour of these turning points when (A, e) varies. 
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4. Reduction to a finite dimensional analysis. Neural field models are one of the pos- 
sible generalizations of standard neural networks considered as discrete sets of connected 
neurons. They can be characterized by two limit processes. First, we let the total number 
of neurons grow to infinity so that each node of the network represents an ideally infinite 
number of neurons, in practice a large number of such neurons belonging to different popula- 
tions, in effect a neural mass. Second, we assume that these neural masses form a continuous 
neural material and let the connectivity graph of the neural network become continuous. The 
graph connectivity matrix then turns into a continuous function of the spatial coordinates. 
One would think that after passing twice to the limit the resulting system would be infinite- 
dimensional. However, the dimensionality of the neural field models depends essentially 
upon the linear operator J representing the connectivity function. If this linear operator has a 
finite-dimensional range, we show below that the corresponding neural field model is finite- 
dimensional and is equivalent to a finite set of ordinary differential equations. Moreover, 
even if this condition is not met we also show that we can always approximate the operator as 
accurately as desired by a finite-dimensional range operator and reduce the neural field model 
to a finite set of ordinary differential equations. 

4.1. The Pincherle-Goursat Kernels. In our numerical studies we use connectivity 
functions that are such that the corresponding linear operators of T have finite rank, i.e. their 
range is a finite dimensional subspace of T. This is without loss of generality because of the 
following theorem (see, e.g. ifTTl '): 

THEOREM 4. 1 . The subspace 1Z f (J-) of finite-dimensional range linear operators of T 
is dense in Ti., the set of linear compact operators of T. 

Proof. This is true because T is a Hilbert space [11, Chapter 6] □ 
In the area of integral equations, these operators are called Pincherle-Goursat kernels PT1 . in 
short PG-kernels. They are defined as follows 

DEFINITION 4.2 (Pincherle-Goursat Kernels). The connectivity kernel J(r, r') is a PG- 
kernel if 

N 

J(r,r') = ^X fc (r)<g>y fc (r') 
fc=i 

where X k , Y k , k = 1 ■ • • N are two sets of N linearly independent elements of T, and 
Xfc(r) ® Yk(r') is the rank 1 p x p matrix X k (v)Y^ (r'). 
We have 

k 

where (•, is a short notation for the inner-product in L 2 (f2, W). Thus J-Ue Span(Xi, ■ • • , 
which we denote by R(J). 

4.2. Persistent state equation for PG-kernels. We now cast the problem of the com- 
putation of the solutions of equation ( |2.10| ) into the PG-kernel framework: 

V - I ext = J • S(AV) 

Since V — I cxt S R(3), we can write V — I cxt = v k X k , and note v = (vk)k=i---N- The 

k 

persistent state equation reads: 



v k = ( Y k , S A [J2 v k X k + I ext ) k = !,■■■ ,N 
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This is a set of N nonlinear equations in the N unknowns v\ , • • • , vn which can be solved 
numerically using classical methods. 

4.3. Reduction to a finite number of ordinary differentiai equations. In this section 
we reduce equation ( |2.1[ ) (when J is independant of t) to a system of ODEs. We write I 
instead of I oxt for simplicity and note Sa the function K p — » MP defined by Sa(x) = S(Ax) 
for x eW. 

We note i?(J) ± the orthogonal complement of R(J) in T with respect to the inner- 
product (, ) 2 : 

T = R(J) © R(J) 1 - 

We write 

V = v" + v- 1 , 

where V" (resp. V- 1 ) is the orthogonal projection of V on R(J) (resp. i?(J)- L ). We have a 
similar decomposition for the external current I 

I = I 11 + I- 1 

We now decompose i?(J) as a Cartesian product of p finite dimensional subspaces of Q. 
Because 

N 

J ij (r,r') = '£xi(r)Yi(r') i, j = 1, ■ ■ • ,p, 

k=l 

each coordinate Vi,i= 1, ■ • • ,p of V satisfies 

JV 

V< + = ^ (n, S(AV)) 2 + U i = 1,. • • ,jj 
fc=i 

Let us consider the p finite dimensional subspaces i = 1, • • ■ ,p of G, where each Ei is 
generated by the N elements XI , k = 1, ■ ■ ■ , N. We note Ef- the orthogonal complement 
of Ei in Q. This induces a decomposition of T as the direct sum of the cartesian product 
Yl? =1 Ei = R(3) and its orthogonal complement Yli=i E t = ^(J)^- We write Vi = 
V)} + V- 1 as well as E = if + I±. We then have 

Vl + ai v! = E <F fe ,S(AV)) 2 + /| 

fc=i »=!,••• ,p (4.1) 

Considering the canonical basis gj, i = 1, • • • ,p, of W, we define 

»=1 i=l 
t=l z=l 

We obtain the 2p-dimensional non-autonomous system of ODEs: 



I Vll + L-V" = J-S(AV)+F 
jv-L + L- V 1 - = I ± 

Remark 9. If I is stationary then converges to L -1 !- 1 . 



(4.2) 
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5. One population of orientation tuned neurons: the ring model. As an application 
of the previous results, we study the ring model of orientation tuning introduced by Hansel 
and Sompolinski (see E7J SDJ \M \TJ\ |7J ED), after the work of Ben-Yishai (see 0), as a model 
of a hypercolumn in primary visual cortex. It can be written as: 



rA(x,t) = -A(x,t) + S | A | / J(x -y)A(y,t))dy/n + I(x) ■ 

tt/2 

Some authors, [7, 9|, chose J to be a difference of Gaussians. On the other hand, Ben- 
Yishai, in 0, started with a network of excitatory /inhibitory spiking neurons and derived 
a meanfield approximation of this network yielding the activity response described by the 
following equations: 

x/2 V 

tA(xA) = -A(x,t) + S ( A ( / [J + Ji cos(a(y - x))] A(y,t))dy/n + el(x) -0 

J(x) = 1 — + /3cos(a(x — xo)) 

a = 2, < /3 < 1 and the threshold 9 = 1 in the above cited papers. Being an activity 
model and not a voltage model in the terminology of lTT8l |2D it is not directly amenable to 
our analysis. We can either extend this analysis to activity models as shown in appendix [Djor 
do the following. We rewrite the previous equation as 

tA = -A+ S(X(J -A + eI-6)), 

and perform the change of variable V = J ■ A + el — 9. This leads to the following equations 

tt/2 

rV(x,t)=-V(x,t)+ J [J + J 1 cos(a(y-x))]S(\V(y,t))dy/ir + eI(x) + 9 

-tt/2 

^J(x) = 1 — (3 + /3cos(a(x — x )) 

(5.1) 

We are now in the case of the model studied in this paper with p = 1, d = 1 and = 
(—tt/2, tt/2). Note that since the x-coordinate represents an angle, this is not a neural field 
model per se but rather a neural mass model. 

The nonlinearity is often chosen to be a Heaviside function, or, as in 0, a piecewise 
linear approximation of the sigmoid]^] or, as in ifTHl [TTl 171. a true sigmoidal function. J% 
can take any values and / is an external current coming from the LGN. Jo is most of the time 
negative (see J31[l7j 7,9]) but can be positive as well(see 13): the Js can be thought of as the 
first Fourier coefficients of J, J , being its mean value, it can be positive even if the surround 
is inhibitory. We can, up to a multiplication of the previous equation, make the assumption 

Jo = e e {-1,1} 

For example, in ifTTl . we find Jo = —7.3, J\ = 11, j3 = 0.1, 9 — which are taken 
from [4] except for 9 — 1. The slope is assumed to be A = 1. Using the previous scaling, it 
becomes Jo = — 1, J\ = 1.5, A = 7.3/si = 29.2 and 9 — > 9/7.3 which gives 9 sa 0.1 in 
the case of and 9 = in IfTTl . 



7 This does not allow for computation of bifurcation branches but allows the detection of branching points. 
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The goal of this section is not to derive the whole bifurcation diagram of the Ring Model 
but rather to show how the stationary solutions are organized and to give clues about the 
dynamics in a given range of parameters. This study is helpful because many large scale 
model of VI use the Ring Model for the hypercolumns. We will see that, depending on 
the stiffness of the nonlinearity, there may exist many stationary solutions, which are all 
acceptable responses of the network for a given input of the LGN. Thus these local orientation 
detectors may behave less trivially than they were initially made for. These cortical states can 
make the local dynamics sophisticated when the slope parameter A is big enough. 

5.1. Mapping the ring model to the PG-kernel formalism. Expanding the cosine in 
the previous equation, and denoting by cos a (respectively sin Q ) the function x — » cos(aic) 
(respectively x — > sin(ax)), we find that, depending on the sign of Jj, (e, = sign(Ji),i — 
0,1): 

2 

J = £ol®l+El \f\-h\ COS a Ji| COS Q +E\ y/\ J\ | Sm Q Jt\ Sm a = ^ SjXjiglXi, S 2 : 

i=0 

This formulation has the advantage of preserving the symmetries of W. With the notations 
of the previous section, we have I 1 - = 0, and 



-V 11 + W ■ 5(A(V r|1 + Vtf-e - */ 7 ")) + e/ 11 + i 



where 



V^(x,t) = vi(t) + v 2 (t)^\J\\cos a x + v s (t)\/\J\\sm a x, 



(5.2) 



i.e. the model is three-dimensional. Note that the previous equation is a rewriting of (5.1 
without any approximation. 
Similarly we have 




(5.3) 



As V 1 - (t) — > 0, we restrict the study to the case V 1 - — even if we lose some of the 'real' 
dynamics by doing so. This is motivated by the fact that the dynamics is made of heteroclinic 
orbits (as we will see in a moment) between persistent states belonging to the vectpr space 
V 1 - = 0. Hence, using this simplification, we are led to study the following 3D system: 



tvi = -V! + £ (S(XV), 1) + ell + 9 
tv 2 = -v 2 + si VW (S(XV) , cos a ) + el\ 
L tv 3 = -v 3 + ei y/\JT\ (S(XV) , sm a ) + ell 



(5.4) 



t/2 

where (/, g) = J f(x)g(x)—, V is given by equation \5.2\ and /]', i = 1, 2, 3 is given by 

-7T/2 

equation ( |5.3[ ). Note that the basis (Xq, Xi,X 2 ) is not orthogonal for this inner product. 

This system enjoys the symmetries described by the following lemma. 

Lemma 5.1. When j| = 0, if\ = {v\ v 2 Vs) is a solution, then so is (v\ v 2 — v 3 ). The 
plane w 3 = is invariant by the dynamics. 

Proof. This is a consequence of the fact that sin a is an odd function while cos a is an 
even function. □ 
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It is easy to see that 

llvll 2 1 f n/2 - dx I \ 

E(v) = -— 1 — - / So(Xvi£o+Xv 2 ei\/\Jx\cos a x+Xv3ei^J\Ji\sin a x) \-{sP+0,v i 

2t t\ J_ 7t/2 tt \ / 



where So is a primitive of So, is an energy function for the dynamics, i.e. v = V-E(v). 
Consequently, even for / non spatially homogenous, there are no non-constant periodic tra- 
jectories nor homoclinic orbit^J Moreover, all bounded trajectories are stationary solutions or 
trajectories converging to stationary solutions. Having proven that all trajectories are bounded 
for the neural field equations in proposition 2.6 we have characterized the dynamics. It re- 
mains to compute the stationary solutions and their attraction basins. 

N-l 

Remark 10. We can generalize these facts to PG-kernels of the type J = ^2 £ iXi ® Xi 

i=0 

by choosing E(y) = + { J2 In S (E Ae^X^r)) dr + (I, v) 

k=l i 

5.2. Finding the persistent states. In order to characterize the set B of stationary solu- 
tions, we apply the scheme of section [XT] Hence we study the following equation 

V = J ■ S (XV) + e/ 11 + n(6 + J ■ 5(0)) 

The nonlinearity is the odd function : 

Sq(x) — z l\ 

Note that J ■ 5(0) = | (e + Ji 2ain °^ /2) cos a ) = \e Q + e x j\J\\ sin ± /2) X^. This gives : 



tv x = -vi + e (So(W), 1) + ell + fx(6 + f ) 

ri 2 = -v 2 + eiv^Jij (5 (AT/),cos a ) + el\ + neiy/W V^"^ (5.5) 
rv 3 = -v 3 + e 1 ^\jT\(S (XV),sm a ) + £ll 

5.2.1. The simpler case [i — e — 0. This corresponds to finding the persistent states 
when (i = e = ensuring that v = is a solution. For the sake of simplicity, we reduce the 
study to the case which breaks the translation symmetry so that we do not get involved 

with equivariant bifurcation theory. 

The Jacobian at v = is given by (using some symmetries): 

-Id 3 x3 + AsiK, 
where sj = 5q (0) = | and the matrix K is equal to 




Ji(l,sin2) _ 

K has in general (for a « 2) three real eigenvalues. Indeed the operator J is self-adjoint 
for the inner-product defined above, but, as previously mentioned, the basis (Xq, Xi,X 2 ) is 



K 




| Ji| (l,cos a ) 
Ji (Mos 2 ,) 




8 This follows from the consideration of the time derivative of the energy E. 
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3 



Table 5.1 

Number of bifurcated Pitchfork branches from (0, A) depending on the values of £q, ei. The value of a in 
\5.\) is close to 2. 



not orthogonal for this inner product, and hence the matrix K is not symmetric in this basis. 
We note <j\ the eigenvalue of K corresponding to the eigenvector (0, 0, 1), a 2 and 03 the 
two eigenvalues ot its upper lefthand 2x2 submatrix. The values, noted A^, i = 1,2,3, 
corresponding to potential n bifurcations are equal to 4/er,. The signs of the A^s give the 
number of bifurcated branches (recall that A > 0). Because s 2 — Sk_(0) = and S3 = 
(0) = — 1/8 7^ 0, all branches are Pitchfork branches (see section |3.l[ ) whose third order 
term is X3 = ( e h e i) 2 (~ H^W e i Hi < for a ~ 2 )- Hence these branches are 



directed toward A > A^. This is summarized in table 5.1 The eigenvectors a, i — 1,2, 3, of 
the Jacobian of j5A\ at v = ar given by: 

ei = sin Q , e 2 , 3 = a 2 , 3 + b 2 ,3 cos a 

We reduce the number of possibilities by assuming from now on that J\ > 0. It turns out 
that in this case a 2 < and there are only two possibilities to consider: o\ < 03 for a > 2 
and 03 < o\ for a < 2. This gives the relative position of the different bifurcated branches, 
noted Pi, i = 1,3. Once we have found the bifurcation point, we can numerically compute 
the bifurcated branches for all positive values of A using a continuation method (we used the 
pseudo-arc length method as described for example in ||39ll32l ). 

From our numerical experiments we conjecture that in the case Eo = —1 and E\ = 1, P\ 
and P3 satisfy the following properties 

1. Pi lies on the t'3-axis. 

2. P3 lies in the plane of equation v$ = 0. 

3. Pi and P3 do not intersect. 

Remark 11. We can reverse the orientation of the pitchforks by choosing a nonlinearity 
such that s 3 > 0, the bifurcation diagram would be more complex: more saddle node points 
would appear because of proposition \3.3\ It is the fact that s 2 = Ofor the sigmoid which pro- 
duces pichfork branches. Another choice of nonlinearity, for example S(x) — i +exp ^_ x+e y 
would produce transcritical branches. It is indeed difficult to imagine that s 2 = for an 
experimental fit of a "real" rate function. But anyway we are not yet looking at the "real" 
bifurcation diagram since we are assuming e = p, = 0. 

Figure |5.1| shows a typical example corresponding to the values of the parameters that 
are found in the largest number of published articles. We come back to this choice in section 



5.3 The left part of the figure shows the three components of the persistent states as functions 
of A. For Pi there is only one nonzero component, U3, in blue. For P3 there are two nonzero 
components, t>i shown in red and v 2 shown in green. The right part of the figure shows 
another representation of Pi and P3 as curves parametrized by A in the (vi,v 2 , U3) space. P3 
is clearly in the (vi,v 2 ) plane while Pi is along the «3-axis. The color at each point of the 
curves represents the value of A according to the color scale shown on the right. 
In detail we have 



'The upcoming nonlinear analysis will show that they are indeed bifurcation points. 
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< A < Ai When A goes through A3, the 0-solution loses its stability and becomes a sad- 
dle. There are three persistent states, (unstable node) and two located on the pich- 
fork branch P3, both stable. The corresponding dynamics is shown in the left part 



of figure 5.1 



Ai < A When A goes through Ai, the 0-solution loses its stability along the «3-axis. There 
are two new persistent states located on the pitchfork branch Pi, both are unstable 
nodes (the unstable manifold is one-dimensional). The corresponding dynamics is 



also shown in the left part of figure 5. 1 
There are at most 5 stationary solutions. 




.0 4 i— . 1 — — 1 . 1 1 1 . 1 — 

2 3 4 5 6 7 8 9 10 



FIG. 5.1. Left: Plot of the bifurcation diagram for a = 2.2. It shows the two pitchfork branches P\ and 
P3. For each branch, we have only plotted the nonzero coordinates with v\ = red, V2 = green for P3 and 
113 = blue for Pi. We have also plotted the dynamics in two and three dimensions according to the values of the 
slope parameter A. Right: Plot of the equilibrium points. The color encodes the value of the slope A (see text). 
J = -1, J 1 = 1.5, u = 0, e = 0, a = 2.2 



Remark 12. In the case Jq = 1 and £\ — 1, there is another pichfork branch to handle, 



see table 5.1 For A big enough, v = becomes an unstable node and there are 7 stationary 



solutions instead of 5 in the case Jq = — 1 

5.2.2. The case \i = 1, e ^ 0. We are now halfway from our scheme completion. To 
have an idea of the persistent states at low contrast (ie e w 0), we need to know the persistent 
states for : A, fj, = 1, e = 0, that is we need to know the solutions of 

vl = J-S(\vl)+6 



Following our program, we numerically compute the persistent states when the slope A 
and fi both vary. As described in section 2.3 we expect many of the previous bifurcations 
to disappear thereby breaking some of the connectivity of the sets B\ which can actually be 
(partially) recovered by considering the sets This was done using the library TRILINOS 
(see l39l and the website ) using multiparameters continuation. 



We show an example of this continuation in figure 5.2 Left where we display the V2 
component of the persistent states as a function (sometimes multivalued) of A and fi. A cross- 
section of this set by the plane of equation /i = (shown as semi-transparent in the figure) 
yields a curve identical to the one shown in green figure 5.1 Left. The figure nicely shows 
how the first pichfork bifurcation branch P3 opens up when ji becomes non zero: this gives 
the connected component of which is linearly stable. 



24 



R. VELTZ AND O. FAUGERAS 




^3 



\ 



FIG. 5.2. Le/f: The V2 component of the 2-parameters continuation (\, p) for 6 = 0.1, o = 2.2, Jo 
— 1, Ji = 1.5. Right: Plot of the equilibrium points. The color encodes the value of the slope X. Jq = — 1, Ji 
1.5, a* = 1, e = 0, a = 2.2, 6 = 0.1. 



O = Ui-£o<So(AF),1) + m(0-, , 

= v 2 - ei VPil (So(AVO,cos Q ) + 
= « 3 -e 1 ^J^(S'o(AF),sin ct > 



»W2) 



(5.6) 



However, as can be seen from |5.2| Right, non-zero values of /i do not break the pitchfork 
Pi . It is easy to qualitatively understand why, even though a full mathematical proof is hard 
to come up with: /i does not affect the third equation which produces the pichfork i\, We 
can prove it locally for /i near using the implicit function theorem. We are looking for a 
point V2(ft), 0) at which a pichfork occurs for A = Ai(/x). Let us consider 



H{v u v 2 , A; n) 



v 1 -e (S (XV),l)+ft(e- 

V2 ~ £l\/JT(5o(AV r ),COS Q ) + fJ£iy/\Ji\ 

1 - AJ X (DSoiXvxXo + Xv 2 X 1 ),sin 2 a ) 



2 > 



«(V2) 



(5.7) 



where the last component of H(v\, v 2l A; zx) is g^(5.6 3). It is easy to see that H (0, 0, Ai, 0) 



[0 0]. The Jacobian of H w.r.t. (wi,w 2 , A) at (0,0, Ai,0) is (because 5^(0) = 0) : 



-1 + AiS^o 

AisieoVPil (1j c °So 




AiSieiVpaJ (l,cos Q ) 
-1 + AisiJi (Lcos 2 ,) 




-s x Ji (l,sin 2 ) 



G GL 3 



Hence there exists a unique solution defined locally for fi > satisfying H(vi(fi), v 2 (ft), Ai(/x); /x) 
: we have found a bifurcated point. Moreover, as Xi (A* = 0) 7^ 0' it w iU remains so for 
small zx: the Pichfork Pi is not affected by /i. For large values of /i we have to rely on 
numerical simulations. 



Now, because of lemma 5.1 the solutions (y%, v 2 ) corresponding to v% (ie lying on 
the pichfork branch) are the same for V3 and — V3 which gives the branch 2-3 in figure [B3| 
Hence, when /1 7^ 0, we still have the pichfork Pi and the branch (located in u 3 = 0) arising 



from the opening of P3. This gives the diagram shown in figure 5.3 for the three components 
of v. 
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-0.4 



' 4 5 6 7 8 9 10 11 12 13 -0.4' 4 _ _l ^ g 

FIG. 5.3. Plot of the three components (v\ ,1)2,113) = (red, green, blue) as functions of the slope Xfor the 
following values of the parameters: Jo = —1, Ji = 1.5, /U = 1, e = 0, a = 2.2, 6 = 0.1. Recall that the slope 
of the model has to be around A = 29. 



This diagram is a bit misleading because if we count the red components, there are 
four of them for A > Ai which gives an even number of solutions (in contradiction with 
proposition|2.9[). In fact by doing so we miss the symmetry v$ <-> — v$ and the corresponding 



3.1 



solutions. It is eas ier to look at figure [572] Right to count the stationary solutions 
From section 



it follows that the connected branch of Vj is stable as well as the 
branch Pi . The only unstable branch comes from the opening of P 3 and is shown in figure 



5.2 Right. 



Remark 13. Figure 5.2 gives a good example where our scheme allows to detect another 



branch of solutions which is not connected to the connected component of the trivial solution. 



The figure 5.3 also tells us which branches will appear when the contrast satisfies e ^ 0: 



this will be a perturbation of figure 5.2 Moreover, our scheme detects two saddle-node at 
A = 20 and A = 26 which are undetectable by A— continuation only. Each saddle-node gives 
2 additional persistent states : one stable, the other unstable. Except in the case xq = 
(see equation ( |5.1| l), the pichfork Pi will open up, giving two other connected components 
in a ddition to the connected component of Vq. As an example, we have plotted in figure 



5.4 



all the persistent states (in the activity representation, ie A? = S(V^), in the case / 
1 — 0.1 + 0.1 + cos a and e — 0.1, A = 9. There are five of them as predicted from a 
perturbation of figure |5.3| 



5.3. Discussion. There are two reasons why we presented this example. First, it is a 
nice simple model to which the formalism of this article easily applies and allows us to push 
the analysis far enough to grasp an almost complete understanding of its persistent states 
and a somewhat detailed understanding of its dynamics. Second, it conveys information for 
models of VI that is likely to be biologically relevant. For example, as the slope A of the 
sigmoid is increased, many (up to 4 depending on the signs of ( Jo, Ji)) new stationary states 
appear whose stability evolves with A. One of these solutions is "dramatic" for the purpose 
of orientation detection: even if the LGN input orientation peaks around the angle xq = 0, 
the Ring Model can produce a stable cortical state (or a percept) corresponding to an angle of 
tt/2! 

However these solutions may be destabilized by adding lateral spatial connections in a 
spatially organized network of ring models; it remains an area of future investigation. As far 
as we know, only Bressloff and co-workers looked at this problem (see (HIED): they studied the 
local dynamics around the bifurcation points but did not look at non-local dynamics (basically 
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they considered that there were 2 stationary solutions around a bifurcation point and forgot 
about the third one predicted by proposition |2.9| l. 

There is no biologically motivated restriction on the values of the slope A which can be as 
large as desired, without mentioning the fact that neural mass models are very often written 
with a Heaviside function for the nonlinearity which, as mentioned previously, is the limit 
case: H(x) = lim S(Xx). 

A — >oo 

We also made the assumption a ^ 2 in the previous analysis. It remains to know how 
much of what preceed holds in the case a = 2. More generally, what would remain if 
one were to choose a difference of Gaussians as a connectivity function over a cortex = 
(— 7r/2, 7r/2). We expect a lot of similarities if the width of the difference of Gaussians is of 
the same order of cos a . 

6. Two populations of spatially organized neurons. We now apply the previous anal- 
ysis to a system we started to analyze in ll23ll . In the somewhat reduced form we consider 
here, it consists of two populations (p = 2), one excitatory, one inhibitory, distributed over a 
fiat (d = 2) cortex — [—1, l] 2 . The connectivity matrix kernel writes: 



J(r,r') = 



aGu(r-r') 
6G 21 (r-r') 



-6Gi 2 (r-r' 
-cG 22 (r-r' 



(6.1) 
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where the Gij(r) = e are two-dimensional Gaussian functions defined on M 2 with 
G12 = G21. a, b, c > characterize the strength of the connections. We also assume 
lext = 0. The parameter [i controlling the translation of the sigmoid (see ( |3.3| l) is therefore 
the only parameter, outside A, that we will vary from to 1. We also chose (notice that 

s 2 + 0) : 

In 11231 . we were able to compute the stationary solutions V{ x when the slope A 
was small (i.e. A < A*) using the Nystrom method. We now know that if we perform a 
continuation of these solutions with respect to A, we are bound to miss quite a few of them. 
Therefore we perform a two-parameter continuation with respect to the pair (A, ji) in order to 
recover more, if not all, stationary solutions. 

Biologically speaking, no systematic investigation has been performed to test the valid- 
ity of the translation invariance of the connectivity function. Hence a roughly translation 
invariant (called heterogeneous in |19|) is not less biologically relevant. This is where the 
PG-kernels are useful: they provide an easy way to approximate the convolution operation as 
well as an effective representation of the connectivity (see section|4|i. 

There are four reasons why we think this example is interesting: 

• We want to show how to deal with heterogeneous kernels. 

• We want to give a non trivial example to the existence of a branch of solutions not 
connected to the trivial solution Vq . 

• We want to show an example of application of proposition 3.3 

• More importlantly, we want to show how the results of section 3.1.1 may change in 
the two-dimensional case. 

As in the previous case, we are not interested in having the complete bifurcation diagram 
of the system but rather in giving numerical examples of the previously enumerated points. 

6.1. PG-kemel approximation of J. Let us write 



V* = e - 



||r||V2 e -||r'|| 2 /2 e <ry> ^ g _ || r f ,2^ || r> f ,2 ^ + ^ + l( r ,r'} 2 + ...) 



We notice two important facts 

1 

2 



• 1 + (r, r') + i(r, r') 2 + ... is a polynomial in the components of r and r', hence a 



PG-kernel. 

s e Hl r ll I 2 has a bell-shape which we approximate with the following function ^1 — 

)a 2 
with a > 0. This choice is also motivated by the fact that e~H r H / 2 tends 

to zero at the edges of the infinite cortex that is usually considered in the litera- 
ture. We keep this property with our finite cortex since ^1 — ||r|| 2 ^ = on the 

boundary d ([— 1, l] 2 ) of [— 1. 1] 2 . Any other positive bell shaped function would be 
appropriate. 

Putting these two facts together, we end-up with the following approximation of a Gaussian 
convolution kernel with a PG-kemel. 

J y (r,r') = dj + (l - ||r|| 2 ) Q, '(l - ||r'|| 2 ) a %(r,r'), i, j=l,2 (6.2) 
where P is a 2 x 2 matrix with polynomial entries in r and r' and C is a constant 2x2 matrix. 



Recall from section 



2.1 



that for d = 2, we require m = 1, i.e. T = W ' (0, Mr). Thus 



J||jf < oo imposes 2(aij — !)>—! hence a^j > 1/2, similarly for a'^. 
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6.2. Numerical experiments. In these numerical experiment we choose an = a' n = 

3, ai2 = a' 12 = a 21 = a 21 = 2, a 22 = & 22 = 4, a = 10, = 15, c = 12.75 and 

<r.r') 

L = 0.5/d. The polynomials Pjj (r , r') are equal to the 4th-order Taylor expansion of e "« , 
with a' 1 = a y . 

This results in the following expressions for the two components v( and V/ of the 



persistent state 



y/ = Q u (r)(l-||r|| 2 ) 3 + Q 12 (r)(l-||r|| 2 



F/ = Q 21 (r)(l-||rf ) +Q 22 (r)(l-||r 



The polynomials Q%j, i, j = 1, 2 are of degree 4. The total number of parameters needed for 
representing a persistent state is therefore equal to 5 x 5 x 4 = 100 variables. This number 
comes from the fact that we have four polynomials Qj», i,j = 1, 2, of total degree 10 in two 
variables that are the products of two polynomials in one variable of degree 4. The reason for 
this is the special form of the polynomials P. L j (r , r') which are functions of the inner product 
(r, r'). This is a somewhat crude approximation of the convolution kernel (6.1 

iach of these 100 va 
jjL. We represent in figure 



Each of these 100 v ariab les is a function of the slope parameter A and contrast parameter 
the infinity norm of the 100-dimensional vector V-^ as a function 
of the slope parameter A for fj, = 0. 



All our analytical predictions are performed, for convenience and simplicity, on the con- 
volution kernel J. They are likely to carry over to the heterogeneous approximate case. 

The general analysis having been conducted (see section |3TT) , we wish to do some spe- 
cific computations to explain the numerical results : the first thing to do is to find the eigen- 
elements of J. 

For each n = {n\, n 2 ) G N consider the function cos n : M 2 — > R defined by cos n (r) = 
cos(niri + n 2 r 2 ), and the vector subspace X COSn of T defined by X COSn — {V | V = 



cos n v, v S 



'}. Xjiii,, is defined in a similar fashion. 



It is easy to check that the vector subspaces X COSn of T are invariant by the linear operator 
J, and, by parity, that J • X siDn = for all n e N 2 . We have 



J ■ (cOS n v) = COS n J(n)v = COS r 



oGi(n) -bG 2 (n) 
bG 2 (n) -cG 3 (n) 



Vv G 



where J, Gi are the Fourier transforms of J and G^. 
Keeping the notations of sectiorjXT] we find that 



U ^(j( m ) 



meN 2 

and the corresponding eigenvector^jare in X COSn for some n G N 2 , hence of the form cos n v 
for some v G M 2 . We have chosen the Gaussians Gi such that the first eigenvalues are simple, 
hence each point A n is a bifurcation point: the corresponding c/z/-factors are 



( n ) / 2 

X 2 (cos„,cos 



3 n ) = 0ifn^0 



10 We do not prove that we obtain all of them though it appears that these particular eigen-elements are sufficient 
to explain the numerically observed facts. 
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A 




FIG. 6.1. Plot of the infinity norm || ■ ||oo of the 1 00-dimensional vector ~V? as a function of the slope parameter 
Xfor fj, = 0. Notice that = is solution for all values of X. The dotted branch ( red and green) connects A ni to 
An 2 ■ The dashed-line (blue) is a set of solutions that is not connected to the component corresponding to the trivial 
solution \f = 0. SN stands for saddle-node. 



xT A o 

Hence Ao is a transcritical bifurcation point whereas the other points A n are pitchfork bifur- 
cation points. Depending on the vector v, the orientation of the pichfork may change. 



In figure 6.1 (case [i = 0), we see two pichforks at A ni , A„ 2 and the transcritical branch 



at Ao- The first pichfork at A ni (red branch) is oriented toward the decreasing As, thus ac- 



cording to proposition 3.3 a saddle-node (noted SN in the figure) must appear. The same is 



true for the transcritical branch (continuous blue branch). Locally (near Ao), each coordinate 
V/ ,i = 1, 2 is 'flat', ie it looks like cosq. 

The bifurcation points A ni and A„ 2 are connected by the red-dotted branch. It cannot 
be seen from this graph that this is true (because two states with the same norm may be 



different), but it can be checked by looking at the 100 components. We plot in figure 6.2 the 
stationary membrane potential V( of the first population for different values of A along the 
branch connecting A ni to A n , . 

The last interesting fact is that our multi-parameter scheme allows to compute bran ches 
that are not connected to the trivial solution (here \A = 0) as one can see in figure 6.1 



the dashed blue branch cannot be detected by A-continuation of the trivial solution = 
because it does not intersect any branch coming from V* = 0. This, again, cannot be read 



directly from figure 6.2 one has to look at all 100 components. 



Remark 14. We have chosen the parameters a, b, c such that the first three eigenvalues 
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have zero imaginary parts. Numerically, most of the other eigenvalues er n have non-zero 
imaginary part leading to Hopf bifurcations. Other stationary bifurcation may appear for 
large slope values (ie A > 40). 




FIG. 6.2. Plot of the first component VJ of the solution along the red-green dotted branch connecting 
A ni to A„ 2 infigure \6.1\ 

7. Conclusion. In this paper we have pursued the analysis, started in ET1I231 . of a spe- 
cial type of integro-differential equation that appears in neural field and neural mass models 
where we are interested in approximating mesoscopic and macroscopic ensembles of neu- 
rons by continuous descriptions. We call these equations neural equations. Like in these two 
previous papers our approach is based on functional analysis to obtain clean results about 
the well-posedness of these equations from which we can study the questions of existence 
and uniqueness of the solutions. In this article we have enforced more spatial regularity by 
choosing the Sobolev, and Hilbert, space W m ' 2 instead of the space L 2 that we used previ- 
ously. The reason for this change is not only technical, it allows us to study the question of 
the bifurcation of the solutions, while ensuring that the membrane potential (or the activity) 
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remains finite on the cortex. 

In effect these equations depend upon a number of biological or experimental parame- 
ters such as the slope A of the nonlinearity, the connectivity matrix J, or the input I. These 
parameters vary in general in neural populations because of such processes as plasticity and 
learning. It is therefore important to understand how the solutions of these equations vary 
when these parameters change. To this end, we used two theories: the degree theory and 
the bifurcation theory. The degree theory does not require the above spatial regularity be- 
cause it describes the general behaviour of the cortical states as the parameters vary. On the 
other hand, the bifurcation theory describes the precise local behaviours of the cortical states 
as the parameters vary and requires spatial regularity to be able to complete the numerical 
computations. 

We believe that a good model should exhibit bounded membrane potentials which re- 
quires to take a bounded nonlinearity (as opposed to some papers, see for example |4|). The 
degree theory yields the powerful estimate that the number of persistent states has to be odd, 
hence it predicts an additional persistent state in the neighborhood of pichfork/transcritical 
bifurcation points. This extra point is invisible to bifurcation theory, which is a local the- 
ory (This was conjectured numerically in lfl9l ). It may change drastically the dynamics and 
shows that local analysis is not sufficient for the study of the neural field equations. 

We have focused on the description of the stationary solutions of the neural equations 
when varying the slope parameter, and tried to compute numerically the additional persistent 
state given by the degree theory for arbitrary external current and connectivity. As the set of 
stationary solutions may not be connected, we used a multiparameter continuation scheme in 
order to compute non-connected branches of solutions and were able to show examples of 
these (see sections[5]and[6]). Whether these different branches of solutions intersect or are un- 
bounded is still unknown in the general case. However, the scalar case for a one-dimensional 



cortex (p = 1, d = l,I ex t = 0) is almost completely solved (see section 3.1.1 1. This point 
has never been mentioned in the literature to our knowledge. The question of whether we 
have computed all the solutions by using our multiparameter scheme is unfortunately still 
open. 

To sum up, using a bounded nonlinearity for the firing rate introduces new stationary 
solutions that were not predicted before. This suggests that the analysis of neural field models 
(of the visual system or of the memory for example) should be re-examined. This is the focus 
of our current efforts. 
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Appendix A. Well-posedness of operators. We prove proposition 2.3 



Proof. The integral in the righthand side of (|2.2[) exists because for almost all r € fl the 



p 2 elements Wy (r, ■, t), i,j = 1, • • • ,p of J are in L 2 (f2) for all t > and the p coordinates 
ofV(-,t) are in L 2 (fT) for alii > 0. 

For each t > this righthand side defines an element of W m,2 (f2). Because of Fubini's 
theorem it is clear that [3(t) ■ V(i)] is an element of L 2 (f2) for all t > 0. 

Next, for each multi-index a, |a| < m, D a [3(t) ■ V(i)] is defined by 

(D«[J(i)-V(t)],$) L2( n) = (-1) W ([J(*) ■ V(t)],D«$) LHU) , 

for all $ in Cg°(f2, W), the space of infinitely differentiable M^-valued functions defined on 
il with compact support. Applying Fubini's theorem again we obtain that the righthand side 
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of the previous equality is equal to 

([D a J(t)-V(t)],$)v {n) . 

This proves that D a [3(t) ■ V(t)] exists for all t > 0, all multi-indexes a, \a\ < m, and is 
equal to [D a 3(t) ■ V(i)]. To see that it is in L 2 (£!) just apply Fubini's theorem again. □ 

Appendix B. Fixed points theorems. We briefly describe some applications of Brouwer's 
and Leray-Schauder's degree theories. They are central in the proofs of parts 4 and 5 of 
proposition '. 



2.9 



We first recall Schaeffer's theorem and provide a short proof based on the Leray-Schauder 
degree. 

THEOREM B.l (Schaeffer). Let X be a real Banach space. Suppose M : X — > X is a 
compact mapping and 

S = {x e X | 3t e [0, 1] such that x = tM (x)} 

is bounded. Then M has a fixed point. 

Proof. We provide for completeness a short proof based on Leray-Schauder's degree 
theory. Taking r > large enough such that S C B°, we define m(x, t) — x — tM(x) on 
B r x [0, 1]. Then ^ m(dB r x [0, 1]) by construction. According to the homotopy invariance 
of the Leray-Schauder degree 

deg LS (Id - M, B r , 0) = deg LS (Id, B r , 0) = 1 ^ 0, 

thus, according to the Kronecker property of the Leray-Schauder degree, there exists a solu- 
tion to the equation M (x) = x. □ 

We can apply this theorem to prove existence of solutions to equation ( |2.10| i. We consider the 
function F : T — * T defined by: 

F(V, A) = -F(V, A) + V = J • S(AV) + 1, 

where F is defined in equation ( |2.10| >. 

It is known that F is a nonlinear compact operator of T [23|. We can apply Schaef- 
fer's theorem to the function F since it is easy to prove (see [23]) that for all V such that 
F(V, A) = tV the following holds 

||V||^<i(AV5jn|||J|| F + ||I||^), (B.l) 

where || || F is the Frobenius norm. 

Hence for all A there exists such that 

F(V{,A)=V{ 

The following easy consequence of the proof of theorem |B.1| is used in the article. 

COROLLARY B.2. For each A > there exists an open bounded set IA\ containing S 
(defined in theorem B.l ) such that deg LS (Id — F,1A\,0) = 1. 
We next give without proof a theorem due to Leray and Schauder. 

THEOREM B.3 (Leray-Schauder). Let X be a real normed space, J — [a,b] and M : 
X x J — > X be of the form Id + m with m : X x J — » X compact on X x J. Let 



S = {(a;, A) e X x J : M(x, A) = 0} 
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and for each A G J, let 

Y, x = {xeX : (x,X) G £} 

Assume that S a is bounded and that 

deg LS (M(.,a),W,0)^0 

/or some open bounded set hi Z> S a . 

77ien £ contains a connected component C intersecting S a x {a} one/ which either in- 
tersects Y,b x {6} or z'i unbounded. 

Appendix C. Compact operators with simple eigenvalues. 

PROPOSITION C. 1 . For every m G N*, the set of compact operators with m simple first 
eigenvalues is dense in the set of compact operators. 

Proof. The set 1Z f(T) of finite dimensional range linear operators is dense in the set 
fi of the linear compact operators of T ifTTTl . Thus we only need to prove the theorem for 
TZf(T). Let us consider J G TZf(T), and R(J) = Span(e\, efc), k > m. Without loss 
of generality we assume that the first eigenvalue of J has multiplicity two, i.e. (3i = fti. Its 
corresponding Jordan block is then 

Pi e 
0i 

Then if we define J „ = J+^e2<g>e2, wehaveliirin^oo J n = J. The first two eigenvalues 
(3\ and 0\ + \ of J„ are simple, i.e. J n ei = /3iei and 3 n (neei + e 2 ) — (/3i + -)(neei + e2). 
We can do the same for the other eigenvalues, and define an operator J„ with m simple first 
eigenvalues, a finite dimensional rank, and which is arbitrarily close to J. □ 

PROPOSITION C.2. The same proposition holds for operators of the type Id + J where 
J is a compact operator. 



Proof. It follows from that of C. 1 □ 



Appendix D. Reduction of the activity based model to a finite number of ordinary 
differential equations. We consider the equation for the activity-based model: 

A = -L„- A + S A (J- A + I) i>0 

A(-,0) = A o (0 K } 

We recall that L a ^ L (see lfT8ll ) because they do not have the same biologiocal meaning: 
One is related to the synaptic time constant and the other to the cell membrane time constant. 
We let L Q = diag(ai , • • • , a p ) 

We also recall the PG-kernel decomposition of J = ^ Xk <8> Yfc 

k 

Each of the p coordinates A i7 i = 1, • • • ,p of A satisfies 

Ai + a l A l = S A (J • A + I); i=l,---,p 

Similarly to the voltage case, let us consider the p finite dimensional subspaces Fi, i = 
1, • • • ,p of Q, where each Fi is generated by the N elements Y£, k — 1. • • • ,N. We decom- 
pose Q as the direct sum of Fi and its orthogonal complement for (, ) 2 F^, Q — Fi © F^~ 
and write Aj = A\ + Aj^. We note Yu (respectively Ylf) tne projection from Q to Fi (re- 
spectively to ; ) parallel to F^~ (respectively to Fi). This induces a decomposition of T as 
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the direct sum of F — Yii=x Q and F = Yll=i Fi sucn tnat f° r eacn vector A of F we 
can write A = A" + A- 1 . By construction we also have 



J A 



0. 



and therefore 



Ai + aiAi = l\f S A (J- A f +T) t 



i = l, 



■P 



which is a 2p-dimensional non-autonomous system of ODEs: 

(AH + L a -A" =11" S A (J • AH +1) 
\A ± +L a -A ± =n ± S A (J.All+I) ' 

where Jl" A = (JlJ A;)j =1 ... p and JT^ A = (Jlf Aj)j = i r .. )P are the projections of A on 
FandF-K 

The first equation is a p-dimensional autonomous system of ODEs, which can be solved 
before solving the second one. 

D.l. Lemmas. Lemma D. 1 . For all x, A € K we /iave 

(S(Aa?) - S(0)) 2 < S(A 2 a; 2 ) - 5(0) 

Proof. We set X = Aa; and consider two cases. 
X > 1 We have > e~ x2 and therefore S(X) - 1/2 < S"(X 2 ) - 1/2. Moreover, since 

S(X) - 1/2 < 1, - 1/2) 2 < - 1 /2 and we are done. 

< X < 1 We let X = log y, 1 < y < e. We therefore have 

1 7/ — 1 1 7/ lo SJ' — 1 

5(X)-l/2 = ^ — - S(X 2 )-l/2= 

y ' ' 2y+l y ' ' 2y l °sv + l 

We consider the expression {y ~ 1 ) 2 (y log v + l)-2(y+l) 2 (y 1 °s w - 1) and prove it is 
negative. Because y losv < y it is upperbounded by (y - l) 2 (y l ° sv + 1) - 2(y l ° sv + 
l) 2 (y l °sv - 1) which has the sign of (y - l) 2 - 2{y 2 lo sv - 1). The last expression 
is upperbounded by (y - l) 2 - 20 log * / - 1) = (e x - l) 2 - 2{e x '' - 1) which is 
negative for < X < 1. 



PROPOSITION D.2. 77ie solutions of equation ( |2.10| > satisfy the following inequalities 
for all A > 



< VWllJI 



del' 



= B 2 



as well as 



del R 
— -D2, 



V{ - V r J 



< A IIJI 



B l 
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Proof. The first inequality is a straightforward consequence of equation ( |2.10| i, taking 
the JF-norm and using the fact that < S(x) < 1 for all x £ ]R.http://fidji. inria.fr/biblio 
For the second one we write V{ — Vq = J • (S(AV{) — S(0)), take the JF-norm of both 
sides of the equality and use the Cauchy-Schwarz inequality. We find 



A v 



< 



S(AV{) - S(0) 



S(V{) - S(0) 



L 2 (n,Rp) 



. But since Vx G R, — i < S(x) - 5(0) < \, we have 

L 2 (n,RJ») 2 2 

< \ ^/p|n|, which proves the first inequality. 



The third inequality can be obtained as follows. It is easy to see that S(Xx) — 5(0) < 4| 



for all a; € E and all A > 0. This implies that 
inequality yields the third. □ 



S(V{)-S(0) 



< 



rf 



The first 



Appendix E. An equality for the adjoint operator. We prove proposition 3.2 



Proof. We prove the proposition in the case rn = p = 1. By definition of JJ- we have 



(U,J* -V 



IT 



(j-u,v: 



VU, V 



We note X = JJ- • V and rewrite the previous equation as 

<U,X) 2 +(VU,VX) 2 = (J-U,V)^ 
We next rewrite the righthand side 

(J-U,V)^= (J-U,V) 2 +(V (J-U),VV> 2 
The crucial step is to observe that 

(V (J.U),W) 2 = <U,(ViJ)$.W) a , 
where Vi indicates the derivative with respect to the first variable, i.e. 

(V 1 J^-VV(r)=X:|^(r' ) r)|X (rVr ' 



Let us fix V in T. X is the solution of 

(U,X) 2 + (vu,vx> 2 



(U,/) 2 VUef 



Hence X is the unique weak solution in T of the homogeneous Neumann problem 



j -AX + X = / on Q 
\ = on dil 

/ is the element of L 2 (£l) defined by 

/ = j^-v + (v 1 j);-vv+ 

The regularity properties of the solutions of elliptic equations l20l Chapter 6] imply that X is 
in W 2,2 (f2). Let us denote by V the differential operator, defined in W 2,2 (f2) by 



V 



Id 
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We now choose V = ej-. This implies that X = AeJ- and hence = in dfl. Furthermore 
we have 

XV ■ e* T = 3* 2 ■ e* T + (Vi J)a ■ VeJ- 
Green's formula with the condition -4^- = show that 

(ViJ)5 • Ve> = -J* 2 ■ Ae>, 

and therefore we have 

XD-e*r = r 2 -(V-e*r), 



which shows that T> ■ ej- is an eigenvector of J£ associated with the eigenvalue A. We can 

V-e* :F = e* L2 . 



choose e* L 2 and e*^ so that 



From Green's formula 



f Be* 

(U, = (U, V ■ e» 2 - jf -^U da. 



an 

— e *L 2 
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